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ABSTRACT 

The  present  study  investigates  compressibility  effects  on  turbulence  by  direct  numerical 
simulation  of  homogeneous  shear  flow.  A  primary  observation  is  that  the  growth  of  the 
turbulent  kinetic  energy  decreases  with  increasing  turbulent  Mach  number.  The  sinks  pro¬ 
vided  by  compressible  dissipation  and  the  pressure-dilatation,  along  with  reduced  Reynolds 
shear  stress,  are  shown  to  contribute  to  the  reduced  growth  of  kinetic  energy.  Models  are 
proposed  for  these  dilatational  terms  and  verified  by  direct  comparison  with  the  simulations. 
The  differences  between  the  incompressible  and  compressible  fields  are  brought  out  by  the 
examination  of  spectra,  statistical  moments,  and  structure  of  the  rate  of  strain  tensor. 
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1  Introduction 


Homogeneous  turbulent  flow  has  long  been  viewed  as  a  simple  model  of  turbulence  which 
is  sufficiently  realistic  to  have  fundamental  physical  features  that  persist  in  more  complex 
flows.  For  example,  hairpin  vortices  which  have  been  observed  in  boundary  layer  experiments 
and  simulations  of  channel  flow,  have  been  identified  through  direct  numerical  simulation 
(DNS)  of  homogeneous  shear  flow  by  Rogers  and  Moin  (1987).  Similarly,  the  longitudinal 
streaks  observed  in  wall  boundary  layers  have  been  shown  by  the  simulations  of  Lee,  Kim 
and  Moin  (1990)  to  exist  in  homogeneous  shear  flow  when  the  shear  rate  is  high  enough 
to  be  comparable  to  those  encountered  near  the  wall  in  boundary  layers.  Rogallo  (1981), 
who  was  the  first  to  perform  DNS  of  this  flow,  showed  by  direct  comparison  with  physical 
experiments  that  the  simulated  flow  statistics  were  in  good  agreement  with  experimental 
data. 

In  the  past  several  years,  there  has  been  a  strong  activity  in  the  area  of  compressible 
turbulent  simulations.  The  first  fundamental  study  of  compressible  isotropic  turbulence 
was  conducted  by  Passot  and  Pouquet  (1987)  who  demonstrated  the  presence  of  shocks 
in  the  flow.  Their  results  motivated  the  study  by  Erlebacher  et  al.  (1990)  to  classify  the 
different  possible  regimes  in  isotropic  turbulence  attainable  from  arbitrary  initial  conditions. 
They  quantified  the  conditions  for  the  presence  or  absence  of  shocks.  Other  results  on 
decaying  isotropic  turbulence  can  be  found  in  Sarkar  et  al.  (1989),  Lee,  Lele  and  Moin 
(1990),  and  Zang,  Dahlburg,  and  Dahlburg  (1992).  Results  on  forced  isotropic  turbulence 
are  also  available  (Kida  and  Orszag  1990).  Simulations  have  now  moved  on  to  the  next 
level  of  complexity,  namely,  homogeneous  turbulence  subjected  to  a  linear  mean  velocity 
field.  For  example,  Coleman  and  Mansour  (1991)  have  considered  turbulence  subjected  to 
homogeneous  compression. 

Homogeneous  shear  turbulence  has  a  production  term  which  continuously  feeds  energy 
into  the  system,  which  then  cascades  down  to  the  smaller  scales.  This  flow  exhibits  some  of 
the  properties  of  incompressible  turbulent  shear  flows,  with  the  added  presence  of  acoustic 
waves,  density  fluctuations,  and  dilatational  velocity  fields.  Direct  numerical  simulations  of 
such  flows  can  provide  complete  statistics  not  available  from  physical  experiments.  Com¬ 
pressibility  effects  on  the  flow  can  be  obtained  from  the  simulations,  and  perhaps  more 
importantly,  mechanisms  responsible  for  these  effects  can  be  identified,  and  turbulence  mod¬ 
els  for  compressibility-related  phenomena  can  be  devised.  The  compressible  problem  was 
considered  by  Feiereisen  et  al.  (1982)  who  performed  relatively  low  resolution  643  simula¬ 
tions  and  concluded  that  compressibility  effects  are  small.  Recently  Blaisdell,  Mansour  and 
Reynolds  (1990)  have  also  considered  compressible  shear  flow,  and  identified  eddy  shocklets 
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in  the  flow  for  high  enough  turbulence  Reynolds  and  Mach  number. 

Research  on  the  small-scale  properties  of  homogeneous  turbulence  is  primarily  focussed 
on  understanding  the  intermittency  of  high  dissipation  regions.  It  has  been  numerically 
established  that  the  rate  of  strain  tensor  has  a  preferred  shape  in  regions  of  strong  dissipation 
(Ashurst  et  al.  1987b).  Furthermore,  the  vorticity  tends  to  align  with  the  intermediate 
eigenvalue  (whose  most  probable  sign  is  positive).  While  this  alignment  property  of  the 
vorticity  is  a  consequence  of  a  simplification  of  the  Euler  equations  (Vieillefosse  1984),  the 
shape  of  the  tensor  is  still  unexplained.  Alignment  of  other  vectors  have  been  considered  by 
Ashurst  et  al.  (1987a). 

In  this  paper,  we  use  databases  from  the  DNS  of  homogeneous  compressible  turbulence  to 
clarify  the  sources  of  reduced  growth  of  kinetic  energy  with  increased  levels  of  compressibility. 
This  is  then  used  to  propose  a  new  decomposition  of  pressure  into  incompressible  and  com 
pressible  components.  In  turn,  this  leads  to  an  improved  model  for  the  pressure-dilatation 
term.  Finally,  we  analyze  the  rate  of  strain  tei  .or  after  decomposition  into  solenoidal  and 
irrotational  components.  Differences  between  the  two  tensors  are  brought  out  through  the 
use  of  one-dimensional  probability  functions. 

2  Governing  Equations 

The  compressible  Navier-Stokes  equations  are  written  in  a  frame  of  reference  moving  with 
the  mean  flow  ulf  This  transformation,  which  was  introduced  by  Rogallo  (1981)  for  in¬ 
compressible  homogeneous  shear,  removes  the  explicit  dependence  on  «i(x2)  in  the  exact 
equations  for  the  fluctuating  velocity,  thus  allowing  periodic  boundary  conditions  in  the  x2 
direction.  The  relation  between  the  moving  frame  x*  and  the  lab  frame  x,-  is 

Xj  =  *i  —  Stx  2  ,  Xj  =  x2  ,  X3  =  x3 

Here  S  denotes  the  constant  shear  rate  u\.2-  In  the  transformed  frame  xf,  the  compressible 
Navier-Stokes  equations  take  the  form 

dtp  +  {pui),i  -  St{pu2 0,i  =  0  (1) 

dt{pUi)  +  ( pujui  ),j  =  -p,-  +  r'ijj  -  SpUnSu 

+  St(pu2'ui)}  1  +  Sip, iSi2  -  Str'a,  1  (2) 

dtp  +  u/pj  +  7  pu'jj  =  StU'/p,  1  +  iStpu'21  +  $ 

+  (7  -  1MT*  -  2 StT„  +  SH%n ]  (3) 


2 


p  —  pRT 


(4) 

where  $  =  Tijdut/dxj  is  the  dissipation  function,  r,-j  the  viscous  stress,  u/  the  fluctuating 
velocity,  p  the  instantaneous  density,  p  the  pressure,  T  the  temperature,  R  the  gas  constant. 
The  molecular  viscosity  p,  and  thermal  conductivity  k  are  taken  to  be  constant.  In  terms 
like  u\  j  in  the  above  system,  the  comma  denotes  a  derivative  with  respect  to  the  transformed 
coordinates  x*. 

A  Fourier  collocation  method  is  used  for  the  spatial  discretization  of  the  governing  equa¬ 
tions.  A  third  order,  low  storage  Runge-Kutta  scheme  is  used  for  advancing  the  solution  in 
time.  The  nonlinear  terms  are  dealiased  by  discarding  modes  lying  in  the  highest  third  of 
the  wave  number  range  resolved  by  the  grid.  This  procedure  fully  dealiases  the  quadratic 
and  not  the  cubic  non-linearities  in  the  Navier-Stokes  equations. 

3  Statistical  Moments  and  Spectra 

We  have  performed  simulations  for  a  variety  of  initial  conditions  and  obtained  turbulent 
fields  with  Taylor  microscale  Reynolds  numbers  Re\  up  to  41  and  turbulent  Mach  numbers 
Mt  up  to  0.7.  Note  that  Re\  -=  qX/v  where  q  =  and  A  =  q/yco^  ,  while  Mt  =  q/c 
where  c  is  the  mean  speed  of  sound.  All  the  simulations  have  been  performed  with  7  =  1.4, 
and  Prandtl  number  Pr  =  0.7.  The  computational  domain  is  a  cube  with  side  r.  The 
results  discussed  here  were  obtained  with  a  uniform  1283  mesh  overlaying  the  computational 
domain. 

In  the  case  of  homogeneous  shear  flow,  the  mean  velocity  «,•  varies  linearly  in  space  and 
remains  invariant  in  time.  The  mean  density  p  is  uniform  initially  and  does  not  evolve  in 
time.  However,  due  to  the  mean  viscous  dissipation,  the  mean  pressure  p  and  temperature 
T  increase  in  time. 

3.1  Effect  of  Compressibility  on  Kinetic  Energy 

Figs.  2-5  show  results  from  two  selected  cases  which  have  identical  initial  data  for  the  velocity 
and  thermodynamic  variables.  The  initial  data  is  incompressible,  that  is,  the  density  is 
constant  and  the  divergence  V.u  =  0.  The  initial  pressure  is  calculated  from  the  Poisson 
equation  appropriate  for  incompressible  flows,  and  the  temperature  is  obtained  from  the 
ideal  gas  equation  of  state.  The  variable  parameters  for  the  problem  are  the  shear  rate  5, 
viscosity  /x,  and  the  speed  of  sound  c.  The  two  cases  have  identical  values  for  S  and  p.,  but 
different  values  for  the  speed  of  sound  c  leading  to  different  initial  Mach  numbers  Mt, o-  Table 
1  lists  the  initial  parameters  for  these  cases.  The  simulation  of  Case  1  ran  up  to  Si  =  18 
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when  Mi  =  0.46  and  Re \  =  43,  while  Case  2  was  continued  up  to  St  =  21  when  Mt  =  0.67 
and  Re\  =  41. 

The  Favre-averaged  kinetic  energy  K  is  defined  by  K  —  pui"ui"/2p.  Note  that  the  overbar 
over  a  variable  denotes  a  conventional  Reynolds  average,  while  the  overtilde  denotes  a  Favre 
average.  A  single  superscript  '  represents  fluctuations  with  respect  to  the  Reynolds  average, 
while  a  double  superscript  "  signifies  fluctuations  with  respect  to  the  Favre  average.  The 
evolution  of  the  kinetic  energy  K  as  a  function  of  non-dimensional  time  St  is  shown  in  Fig.  1. 
Physical  experiments  and  DNS  of  the  incompressible  case  indicate  that  K,  after  an  initial 
transient,  evolves  as  Kq  exp(A St).  The  feature  of  exponential  growth  is  carried  over  to  the 
compressible  case,  as  evidenced  by  the  approximately  linear  regime  in  the  plot  of  \n(K/Ko) 
in  Fig.  2  for  St  >  7.  Although  the  picture  of  exponential  growth  survives,  the  growth  rate  of 
K  shows  a  significant  decline  with  a  decrease  in  the  speed  of  sound.  Figure  3  shows  that  the 
growth  rate  A  =  d(lnK/Ko)/d(St )  decreases  by  about  30%  when  Mt  increases  from  0.2  to 
0.4.  We  note  that  the  reduction  of  kinetic  energy  growth  rate  with  increasing  compressibility 
is  a  consistent  trend  in  all  our  simulations. 

In  order  to  understand  the  phenomenon  of  reduced  growth  rate  of  kinetic  energy,  we 
consider  the  equation  governing  the  kinetic  energy  of  turbulence  in  homogeneous  shear  which 
is 

jfiK)  -  pP  ~  P*s  ~  P*c  +  P'd'  (5) 

where  V  =  — 5u"«"  is  the  production,  es  =  the  solenoidal  dissipation  rate,  ec  = 

(4/3)i7d'2  the  compressible  dissipation  rate  and  j/d'  the  pressure-dilatation.  The  last  two 
terms  represent  the  explicit  influence  of  the  non-solenoidal  nature  of  the  fluctuating  velocity 
field  in  the  kinetic  energy  budget. 

If  the  notion  of  exponential  growth  in  homogeneous  shear  flow  is  correct,  A  should  asymp¬ 
tote  to  a  constant.  The  relative  importance  of  the  various  terms  on  the  rhs  of  (5)  cam  be 
gauged  by  their  contribution  to  A.  We,  therefore,  normalize  (5)  by  pSK.  The  terms  on  the 
rhs  of  (5)  are  plotted  in  Fig.  4  after  being  thus  normalized.  According  to  Figs.  4a-4b,the 
dilatational  terms  are  negative  and  reduce  A.  By  St  =  20,  the  combined  contribution  of 
the  dilatational  terms  is  approximately  20%  of  A,  implying  that  they  need  to  be  considered 
in  turbulence  modeling.  The  normalized  production,  V/SK ,  is  plotted  in  Fig.  4c.  The 
normalized  proauction  in  Case  2  begins  to  deviate  from  that  in  Case  1  around  St  =  5,  and 
eventually  asymptotes  to  a  value  noticeably  smaller  than  in  Case  1.  Thus,  the  production  is 
reduced  due  to  compressibility.  The  normalized  dissipation  in  Fig.  4d  also  shows  a  decrease 
in  Case  2  relative  to  Case  1,  but  to  a  smaller  degree  than  the  production.  The  reduction  in 
growth  rate  of  K  in  Case  2  relative  to  Case  1  seems  to  be  primarily  due  to  the  dilatational 
terms  during  the  early  phase  of  the  evolution  and,  for  the  later  phase  (St  >  8),  is  related  to 
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a  decreased  level  of  production. 


3.2  Measures  of  Compressibility 

Apart  from  the  differences  in  the  evolution  of  turbulence  statistics  in  the  case  of  compressible 
flow  relative  to  incompressible  flow,  a  question  that  arises  is  how  compressible  is  the  turbu¬ 
lence?  The  turbulent  Mach  number  Mt  =  y/TK/c  and  the  normalized  rms  density  prms/p  are 
measures  of  departure  from  incompressibility,  since  both  these  quantities  are  zero  for  strictly 
incompressible  flow.  Mt  shows  a  monotone  increase  in  Fig.  5a  reaching  a  maximum  of  0.6, 
which  is  larger  than  the  upper  bound  of  Mt  in  free  shear  layers  and  wall  boundary  layers 
encountered  in  aerodynamic  practice.  Figure  5b  shows  that,  after  starting  from  zero  density 
fluctuations,  both  simulations  develop  significant  rms  density  levels  at  St  =  18,  and  most  of 
the  increase  in  rms  density  occurs  rapidly  within  St  <  3  for  Case  2  (with  the  higher  initial 
Mt).  The  velocity  field  is  no  longer  solenoidal  in  compressible  flow,  that  is,  d  =  V-u  ^  0.  For 
homogeneous  flow,  the  Helmholtz  decomposition  gives  a  unique  decomposition  of  the  veloc¬ 
ity  field  into  incompressible  and  compressible  components  by  u  =  u1 4-  uc,  where  V  •  u1  =  0 
and  V  x  uc  =  0.  The  departure  of  the  velocity  field  from  solenoidality  is  another  measure 
of  compressibility,  and  can  be  represented  by  the  ratio  of  dilatational  variance  to  enstrophy 
Xd  =  d!2 1 or  by  the  fraction  of  kinetic  energy  which  is  dilatational  \K  —  Kc/ K.  The  for¬ 
mer  quantity  \d  is  more  general,  because  the  Helmholtz  decomposition  required  for  obtaining 
Kc  is  unique  only  for  homogeneous  flows.  Figure  6  shows  that  after  an  initial  transient,  Xd 
increases  monotonically  for  Case  1,  but  levels  off  at  about  8%  for  Case  2.  The  quantity  xk 
(not  plotted  here)  is  comparable  to  Xd,  leveling  off  at  approximately  6%  for  Case  2.  Since 
the  mean  shear  leads  to  anisotropy  of  the  Reynolds  stress  tensor,  the  dilatational  fraction 
of  each  Reynolds  stress  component  is  obtained.  Except  for  the  transverse  component  u2u2, 
the  dilatational  contribution  to  the  Reynolds  stress  tensor  was  a  few  percent,  comparable  to 
the  magnitude  of  xk ■  The  ratio  of  u2'u2' /u2u2  is  substantially  larger,  as  shown  in  Fig.  7. 
Since  the  dilatational  component  a  relatively  poor  mixer  of  scalar  compared  to  the  vortical 
component,  the  large  ratio  of  u2‘ 'u2‘ ’ /u2u2  implies  that  mixing  is  preferentially  decreased  by 
compressibility  in  the  direction  of  shear. 

3.3  Compressible  and  Solenoidal  Spectra 

The  Fourier  component  of  the  velocity  is  decomposed  into  components  perpendicular  and 
parallel  to  the  wave  number  vector  from  which  the  solenoidal  spectrum  Es(k )  and  com¬ 
pressible  spectrum  Ec(k)  are  calculated.  Figure  8  compares  the  solenoidal  and  compressible 
spectra  at  St=15  (when  Mt  =  0.55)  for  Case  2.  The  compressible  energy  is  small  relative  to 
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the  incompressible  energy  for  all  modes.  The  compressible  spectrum  is  flatter  at  low  wave 
numbers  ( k  <  16)  relative  to  the  incompressible  spectrum.  The  two  spectra  have  similar 
slopes  in  the  intermediate  wave  number  range.  The  solenoidal  spectrum  is  compared  in  Fig.  9 
between  an  incompressible  run  (Mt  =  0)  and  a  compressible  run  at  St  =  17.  The  compress¬ 
ible  run  has  an  initial  Mt  of  0.4  and  the  initial  velocity  and  pressure  fields  are  the  same  as  in 
the  incompressible  run.  ^From  Fig.  9,  it  appears  that  the  shape  of  the  solenoidal  spectrum  is 
not  altered  by  compressibility,  even  though  the  compressible  fluctuations  are  non-negligible 
at  this  time-  pnns/p  =  0.12  and  \k  =  Kc/K  =  0.05.  However,  the  pressure  spectrum  in 
Fig.  10  shows  significant  differences  between  these  two  cases.  In  the  compressible  case,  the 
pressure  spectrum  seems  to  be  relatively  flatter  than  in  the  incompressible  case. 

4  Modeling  the  Dilatational  Terms 

We  showed  in  Fig.  4  that  the  compressible  dissipation  ec  and  the  pressure-dilatation  p'd' 
contribute  significantly  to  the  kinetic  energy  budget  and  therefore  require  modeling.  In 
Sarkar  et  al.  (1991),  we  proposed  a  model  for  the  compressible  dissipation  ec  ==  aiC,Mf  based 
on  an  asymptotic  analysis  and  DNS  of  isotropic  turbulence.  In  the  present  simulations,  after 
starting  from  a  variety  of  initial  conditions,  ec  ~  0.5 e,Mf2,  suggesting  that  c*i  =  0.5. 

Our  direct  numerical  simulations  of  isotropic  turbulence  and  homogeneous  shear  flow 
provided  a  data  base  for  the  pressure-dilatation  and  suggested  a  theoretical  approach  towards 
modeling  it.  The  evolution  of  the  pressure-dilatation  p'd'  for  Case  1  is  depicted  by  the  solid 
curve  in  Fig.  11.  £From  numerical  experiments,  it  was  found  that  the  nominal  time  period 
of  the  oscillations  in  fid'  decreased  approximately  linearly  with  the  speed  of  sound.  This 
suggested  that  one  could  isolate  the  oscillatory  part  of  p'd'  by  decomposing  the  fluctuating 
pressure  p'  into  the  sum  of  an  incompressible  part  p1'  and  a  compressible  part  pc>.  The 
incompressible  pressure  p1'  satisfies 

VV '  =  -2 uijpu'jti  -  2 uijpu'jj  -  p{ui'u/),ij  ~  2 (6) 

while  the  remainder  pc>  is  the  compressible  pressure.  The  rhs  of  (6)  collects  all  terms  that 
depend  on  p  in  the  equation  obtained  by  taking  the  divergence  of  the  momentum  equation. 
In  the  simulations,  since  p'  is  available  from  the  compressible  Navier-Stokes  solution  and 
p1'  is  evaluated  from  (6),  we  obtain  pc'  as  the  difference  p'  —  p1' .  Since  p'  =  p1'  +  pc'  we 
have  p'd'  =  pl'd'  +  pc'd'.  The  oscillations  are  substantial  only  for  pc'd'  (dotted  curve  in 
Fig.  11),  and  furthermore,  the  peaks  and  valleys  in  the  evolution  of  pc'd'  seem  to  be  much 
more  symmetric  around  the  origin  than  those  in  p'd'  (solid  curve  in  Fig.  11).  The  component 
pl'd'  (dash-dotted  curve  in  Fig.  11)  does  not  have  strong  temporal  oscillations,  and  shows 
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a  systematic  decrease  with  time.  In  order  to  gauge  the  relative  importance  of  the  two 
components  pc'd'  and  pl'd'  of  the  pressure-dilatation  in  the  evolution  of  the  turbulent  kinetic 
energy,  we  calculate  the  time  integrals  of  these  components.  The  integrated  contribution  of 
p1' d'  is  about  an  order  of  magnitude  larger  than  that  of  pc'd'  in  Case  2.  Examination  of 
other  DNS  cases  indicates  that,  in  general,  pc'd'  has  a  negligible  contribution  to  the  turbulent 
kinetic  energy  evolution  relative  to  p1'  d’ .  Therefore,  it  seems  that  only  the  component  p1' d' 
of  the  pressure-dilatation  requires  modeling  in  shear  flows. 

In  order  to  model  pv d\  we  consider  the  Poisson  equation  (6)  for  the  incompressible 
pressure.  After  splitting  the  pressure  into  a  rapid  part  pR<  and  a  slow  part  ps> ,  we  obtain  the 
following  exact  expressions  for  the  rapid  pressure-dilatation  and  slow  pressure-dilatation: 

fi^  =  2pum,nJ^Enjdk  (7) 

ps‘d'  =  pj  -  iulujmum)dk  (8) 

A  A 

Here  <fi  denotes  the  complex  conjugate  of  the  Fourier  transform  <f>  and  Enj  represents  the 
spectrum  of  the  Reynolds  stress  tensor  u'nu'}.  Using  scaling  arguments  (see  Sarkar  (1991)  for 
details)  to  simplify  (7)-(8)  we  find  that  pR'd'  depends  on  the  production  V  while  ps'd‘  depends 
on  the  dissipation  ea.  Finally,  we  propose  the  following  model  for  the  pressure-dilatation: 

_  g 

P'd'  =  a2pUijbijq2Mt  +  azpe,M2  +  - pui,ix(Mt)q 2  (9) 

where  =  u'-u"/q 2  —  6{j/ 3  is  the  anisotropy  tensor,  Mt  =  \2Kjc  the  turbulent  Mach 
number,  and  e,  the  solenoidal  dissipation.  To  obtain  the  functional  dependence  x(iV/()  in 
the  last  term  of  Eq.  (9),  one  would  require  data  from  a  flow  with  mean  homogeneous  com¬ 
pression.  In  this  paper,  we  validate  and  calibrate  the  first  two  terms  in  the  model  for  p?d'. 
In  homogeneous  shear,  the  model  becomes 

fid'  =  -Q.ipPMt  +  azpesM f  (10) 

where  V  =  — Su'{u 2.  Because  the  production  V  =  0  in  decaying  isotropic  turbulence,  the 
variation  of  the  incompressible  pressure-dilatation  with  c3  can  be  verified  using  DNS  of 
isotropic  turbulence.  The  ratio  p1' d' / (pcsMf)  is  shown  as  a  function  of  non-dimensional 
time  in  Fig.  12.  The  decaying  isotropic  turbulence  simulations,  Dl,D2  and  D3  start  with 
Mtfl  of  0.6, 0.5,  and  0.4  respectively.  pl'd'/(pesM t2)  reaches  an  equilibrium  value  by  a  time 
of  0.25,  substantiating  the  validity  of  the  second  term  in  (9).  Based  on  the  DNS  value  of 
the  equilibrium  ratio,  the  model  coefficient  a3  in  (9)  is  taken  to  be  0.2.  The  remaining  part 
of  the  model  for  the  pressure-dilatation  is  calibrated  against  simulations  of  homogeneous 
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shear  flow.  Fig.  13  shows  results  from  three  cases  Sl,S2  and  S3  with  Mt> o  of  0.2, 0.3,  and 
0.4  respectively.  Cases  S2  and  S3  have  shear  S=15,  while  Case  Si  has  S=20.  After  St  =  5, 
the  ratio  (p/,d/  —  0.2 pe3M?) / (j>V Mt)  evolves  in  a  similar  fashion  for  the  three  cases  and 
exhibits  a  slight  decrease  with  time.  In  the  period  5  <  St  <  15,  the  three  cases  S1-S3 
encompass  a  range  of  turbulent  Mach  numbers  of  0.21  <  Mt  <  0.65.  The  results  of  Fig.  13 
suggest  that,  with  a  coefficient  oc2  =  0.15,  the  proposed  model  (9)  is  able  to  parametrize  the 
pressure-dilatation  for  0.2  <  Mt  <  0.6. 

5  Structure  of  the  Rate  of  Strain  Tensor 

In  this  section,  we  study  several  properties  of  the  rate  of  strain  tensor  StJ  to  help  identify 
characteristics  of  the  flow  which  are  solely  the  result  of  compressibility  effects,  and  how 
these  characteristics  differ  from  corresponding  properties  in  an  incompressible  turbulent  flow. 
Because  the  flow  is  anisotropic,  it  is  useful  to  consider  characteristics  of  the  flow  both  with 
respect  to  the  laboratory  frame  of  reference,  and  with  *espect  to  a  frame  of  reference  attached 
to,  and  rotating  with  an  individual  fluid  element.  Global  orientation  properties  of  Sij  are 
established  by  studying  the  principal  directions  of  Sij  in  the  laboratory  frame  of  reference. 
(In  isotropic  turbulence,  this  diagnostic  does  not  provide  useful  information.)  On  the  other 
hand,  local  properties  can  be  characterized  by  considering  the  relative  orientation  of  vectors 
such  as  velocity,  vorticity  or  scalar  gradients,  with  respect  to  the  principal  directions  of  Sij. 
We  consider  the  rate  of  strain  tensor,  decomposed  according  to 

Sij  =  S?j  +  S!j  +  ^ui,i6ij  (11) 

where  Sjj  and  S£  are  respectively  constructed  from  the  solenoidal  and  compressible  (irrota- 
tional)  velocities  u;  and  uc.  Both  Sjj  and  Sjj  are  deviatoric.  The  properties  of  these  tensors 
are  expressed  in  terms  of  invariant  quantities,  i.e.  eigenvalues  and  eigenvectors.  Because  Sij 
is  symmetric,  its  eigenvalues  A,,  i  =  1,2,3  are  real,  and  its  eigenvectors  are  in  the  directions 
of  maximal  or  minimal  extension  (depending  on  the  sigi of  the  A,-)  of  the  tensor.  These 
directions  are  given  by  the  tensor’s  three  eigenvectors.  For  future  reference,  the  eigenvalues 
are  ordered  from  smallest  to  largest: 

A,  <  A,  <  A,.  (12) 

The  first  three  invariants  of  the  tensor  Sij.  defined  as  the  coefficients  of  the  polynomial 
characteristic  equation,  are  related  to  small  scale  phenomena,  and  can  be  expressed  in  terms 
of  the  eigenvalues.  These  relationships  are: 

/  =  —  (Aj  +  A2  +  A3) 
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II  —  Ai  A2  +  Ai  A3  +  A2A3 
III  =  -AXA2A3 

If  the  tensor  is  deviatoric,  1  =  0,  —  II  =  SiAi  which  for  incompressible  flow  is  proportional 
to  the  dissipation,  and  in  isotropic  flow,  III  is  proportional  to  the  third  moment  of  the 
velocity  gradient  probability  distribution  function  (pdf). 

We  present  results  from  one  simulation  with  S  =  15,  and  n  =  1/150,  Mt 0  =  0.3,  R\o  =  20 
and  incompressible  initial  data.  The  simulated  flow  is  analyzed  at  St  =  1,3, 5, 7, 9,  which 
are  all  sufficiently  resolved.  More  details  of  this  analysis  can  be  found  in  Erlebacher,  Sarkar 
and  Hussaini  (1991).  Note  that  all  the  pdf’s  that  foilow  are  unconditioned  and  unweighted. 
Sampling  for  the  pdf’s  (probability  density  function)  is  done  on  a  grid  resolution  of  4£x9Cx96 
although  the  simulation  was  performed  on  a  963  mesh.  This  provides  over  400, 000  sample 
points. 

The  pdf  of  Aj/Aj  and  Af/Af  are  shown  in  Figs.  14  and  15.  Both  plots  show  that 
after  an  initial  transient,  the  eigenvalue  .atios  have  single  peaks.  These  are  located  at  - 
0.25  for  the  solenoidal  ratio,  and  at  approximately  -0.375  for  the  irrotational  ratio.  As  a 
consequence,  the  rate  of  strain  tensors  Si  and  S£,  respectively,  have  preferred  strains  in  the 
ratios  (—4  :  1  :  3)  and  (—8  :  3  .  5)  along  the  principal  axes.  Conditioning  of  A^/ Aj  with 
respect  to  e,  sharpens  the  peak  for  higher  values  of  dissipation.  Further  processing  of  the 
irrotational  ratio  distribution  is  underway.  A  check  of  the  pdf  of  the  irrotational  ratio  was 
also  performed  from  a  128‘J  database  (Case  2)  with  800,000  sample  pointc,  and  its  shape 
is  qualitatively  similar  with  the  peak  at  -0.375.  Note  that  the  pdf  of  the  irrotational  ratio 
is  more  broadband  than  that  of  the  solenoidal  ratio.  However,  the  location  of  the  peak  is 
well  defined,  and  constant  in  time.  This  equilibrium  structure  of  (— 4  :  1  :  3)  for  Sii  is 
observed  in  incompressible  flow,  first  by  Ashurst,  Kerstein,  Kerr  &  Gibson  (1987).  They 
considered  the  statistics  of  (Aj)2/(5^5f )  conditioned  on  the  dissipation  and  found  that  the 
preferred  ratio  was  most  dominant  in  regions' of  high  dissipation.  In  a  later  work,  Chen  et 
al.  (1990)  displayed  scattergrams  of  II  versus  III  based  on  incompressible  mixing  layer 
DNS  data,  which  clearly  demonstrated  that  III  is  proportional  to  (— 7/)3^2  in  the  regions  of 
highest  dissipation.  This  is  in  fact  a  statement  about  the  preferred  shape  of  the  Sj  principal 
ellipsoid.  Our  results  indicate  that  compressibility  does  not  substantially  affect  this  preferred 
structure  of  *he  solenoidal  S-j. 

The  distribution  of  A{  and  Af  show  some  distinctive  differences,  as  illustrated  in  Figs. 
16-17.  Both  figures  show  negative  skewness  of  the  distributions,  and  a  flattening  in  time, 
although  the  effect  is  mu  Ji  more  severe  for  Af  On  both  plots,  the  most  probable  value 
foi  the  eigenvalue  shifts  towards  more  negative  values.  In  the  irrotational  case,  the  effect 
is  extremely  strong.  However,  the  most  probable  ratio  of  Af  to  Af  (see  Fig.  15)  remains 
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invariant  in  time.  Figure  17  also  indicates  that  the  amount  of  change  in  the  distribution  of 
Af  decreases  with  time,  but  there  is  no  strong  evidence  that  the  shape  of  the  pdf  is  evolving 
towards  a  steady  state. 

Information  on  structural  differences  in  the  flow  as  they  relate  to  preferred  directions  of 
straining  for  both  the  solenoidal  and  dilatational  components  of  the  flow  are  presented  next. 
After  the  eigenvectors  of  the  rate  of  strain  tensors  are  normalized  to  unity,  we  compute  the 
angle  fly  between  eigenvector  f,-  and  the  unit  vector  in  the  coordinate  direction  xj.  The  pdf’s 
of  |  cos  flf-  1  and  |  cos 6fj\  are  then  computed.  The  cosine  of  fl  is  chosen  instead  of  0  so  that 
the  probability  density  function  based  on  a  Gaussian  distribution  for  the  velocity  derivatives 
is  flat.  In  Figs.  18-19,  we  respectively  plot  the  pdf  of  jcosflf,]  and  |  cos  flpj,  (i  =  1,2,3)  at 
St  —  9  to  illustrate  that  the  solenoidal  and  irrotational  rate  of  strain  tensors  have  different 
preferential  alignments  in  the  laboratory  frame  of  reference.  Note  that  R\  =  33  at  this  time. 
The  direction  of  maximum  compression  of  Sjj  (f/)  tends  to  be  at  approximately  45°  to  both 
the  x  and  y  axes,  with  no  preferred  orientation  with  respect  to  z.  The  alignment  of  f-f  is  also 
in  the  45°  direction,  while  f/  has  no  unambiguous  preferred  direction.  The  irrotational  rate 
of  strain  tensor  exhibits  markedly  different  properties.  Most  of  the  distributions  of  cos  fly- 
have  two  strong  peaks,  thus  indicating  two  preferred  directions  of  alignment.  For  example, 
fp  is  most  often  aligned  along  one  of  the  tin « •  coordinate  axes.  This  is  true  of  the  other 
two  principal  directions,  but  to  a  lesser  degree.  Interestingly  enough,  all  three  pdf’s  of  if 
exhibit  a  valley  in  the  neighborhood  of  45°.  Statistics  based  on  the  database  of  Case  2  (at  a 
higher  shear  and  higher  Reynolds  number  than  the  results  plotted  in  the  preceding  figures) 
show  that  the  pdf’s  of  j  cos  fl(,|  and  I  cos  flp|  are  sharper.  In  particular,  the  tendency  for  f/  to 
align  is  strongest  in  a  direction  at  30°  to  the  x  axis  and  not  45°.  Finally,  we  notice  that  all 
the  pdf’s  of  ]  cos  fly]  tend  towards  a  steady  distribution  with  increasing  St.  The  distribution 
of  ]  cos  Of- 1  reaches  an  equilibriated  state  at  a  much  earlier  stage  than  does  jcosflhj. 

6  Conclusions 


We  have  performed  DNS  of  homogeneous  shear  using  a  1283  grid  to  ascertain  the  influence  of 
compressibility  on  turbulence  statistics  and  structure.  The  primary  result  Is  the  stabilizing 
effect  of  compressibility  on  the  growth  of  kinetic  energy.  We  find  that  the  reduction  in 
growth  rate  is  due  to  the  pressute-dilatation  and  compressible  dissipation  acting  as  sinks  for 
the  kinetic  energy,  and  due  to  the  reduced  level  of  the  Reynolds  shear  stress.  After  starling 
w*th  zero  density  fluctuations,  the  simulations  develop  significant  density  fluctuations  and 
moderate  levels  of  energy  in  the  dilatational  component.  The  dilatational  velocity  component 
is  strongly  anisotropic,  and  consequently  its  contribution  to  the  transverse  rms  velocity  is 
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much  larger  than  to  the  other  rms  turbulent  velocities.  Although,  the  turbulent  Mach  number 
and  rms  density  fluctuations  are  significantly  large,  the  shape  of  the  energy  spectrum  is 
practically  unaffected  by  compressibility.  However,  the  pressure  spectrum  is  much  flatter  for 
the  compressible  case  relative  to  the  incompressible  case.  The  eigenvalues  and  eigenvectors 
of  both  the  solenoidal  and  irrotational  rate  of  strain  tensors  were  examined  for  preferred 
structure.  In  principal  axes,  the  solenoidal  rate  of  strain  tensor  has  a  preferred  shape  of  (- 
4:1:3)  in  accord  with  previous  results  available  for  incompressible  flow.  The  irrotational  rate 
of  strain  tensor  has  a  different  preferred  shape.  The  pressure-dilatation  correlation  needs  to 
be  modeled  for  compressible  flow.  We  consider  a  form  of  the  pressure  equation  obtained  by 
taking  the  divergence  of  the  momentum  equation,  deduce  a  formal  solution  for  the  pressure- 
dilatation,  and  then  simplify  to  obtain  a  model  for  the  pressure-dilatation.  The  model  seems 
to  compare  well  with  DNS  data  in  isotropic  turbulence  and  homogeneous  shear  turbulence. 
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Fig.  4b.  Contribution  of  compressible  dissipation  to  the  growth  rate  of  K, 


Fig.  4c.  Contribution  of  production  to  the  growth  rate  of  K. 


Fig.  5a.  Evolution  of  turbulent  Mach  number. 


Fig.  5b.  Evolution  of  rins  density. 


Fig.  6.  Evolution  of  ratio  of  dilatational  variance  to  enstrophy. 


log(k) 

Fig.  10.  Pressure  spectra  for  compressible  and  incompressible  runs. 


Fig.  11.  Evolution  of  pressure-dilatation  in  Case  1. 
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Fig.  12.  Incompressible  pressure-dilatation  in  decaying  isotropic  turbulence 


Fig.  13.  Incompressible  pressure-dilatation  in  homogeneous  shear  turbulence. 
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